function [ vSignature ] = compIntrisicWaveDescriptor( distanceMtr, vertices, triangles, desciptDist )
%COMPINTRISICWAVEDESCRIPTOR Summary of this function goes here
%   Detailed explanation goes here

[nbV, ~]=size(distanceMtr);
nbIsoCurv=numel(desciptDist);
[nbTri,~]=size(triangles);

isoCurvLength=zeros(nbV,nbIsoCurv);

global gIsocurvEdges;
gIsocurvEdges = zeros(nbIsoCurv, nbV, 0);


for i=1:nbTri
    isoCurvLength=compFor1Tri(distanceMtr, vertices, triangles(i,:),desciptDist,isoCurvLength);
end

vSignature=isoCurvLength;

% Compute the normalisation: divide by 2.PI.radius
for i=1:size(vSignature,1)
    for j=1:size(vSignature,2)
           vSignature(i,j)=vSignature(i,j)/(2*pi*desciptDist(j));
    end
end

%exportIsocurvEdges(2);

end % compIntrisicWaveDescriptor


function [ isoCurvLength ] = compFor1Tri( distanceMtr, vertices, tri, desciptDist, isoCurvLength)

global gIsocurvEdges;

[nbV, ~]=size(distanceMtr);
nbIsoCurv=numel(desciptDist);
for i=1:nbV
    % Skip triangles that are far away
    if (distanceMtr(tri(1),i)>=desciptDist(nbIsoCurv) && ...
            distanceMtr(tri(2),i)>=desciptDist(nbIsoCurv) && ...
            distanceMtr(tri(3),i)>=desciptDist(nbIsoCurv))
        continue;
    end
    % For each isocurve
    for j=1:nbIsoCurv                
        % Check if the triangle overlap the isocurve
        vInside=[];     % Vertices inside the isocurve
        vOutside=[];    % Vertices outside the isocurve
        for k=1:3
            if distanceMtr(tri(k),i)<desciptDist(j)
                vInside=[vInside,k];
            else
                vOutside=[vOutside,k];
            end
        end
        if numel(vInside)==0
            continue;
        end
        if numel(vOutside)==0
            continue;
        end
        % The triangle overlap the isocurve
        % Compute the two edges which overlap the isocurve
        if numel(vInside)==1
            v1=vInside(1);
        else
            v1=vOutside(1);
        end
        switch v1
            case 1
                v2=tri(2);
                v3=tri(3);
            case 2
                v2=tri(1);
                v3=tri(3);
            case 3
                v2=tri(1);
                v3=tri(2);               
        end
        v1=tri(v1);
        
        % Compute the position of the isocurve along the edge
        ratio=(desciptDist(j)-distanceMtr(i,v2))/(distanceMtr(i,v1)-distanceMtr(i,v2));
        assert(ratio>=0 && ratio<=1);
        p2=ratio*vertices(v1,:)+(1-ratio)*vertices(v2,:);
        
        ratio=(desciptDist(j)-distanceMtr(i,v3))/(distanceMtr(i,v1)-distanceMtr(i,v3));
        assert(ratio>=0 && ratio<=1);
        p3=ratio*vertices(v1,:)+(1-ratio)*vertices(v3,:);
        
        % Compute the corde length
        chordL=norm(p2-p3);
        % Compute the arc length, does not improve the results, make them
        % even worse
        % arcL=2*desciptDist(j)*asin((chordL/2)/desciptDist(j));
        arcL=chordL;
        
       
       curIsoEdges = gIsocurvEdges(j, i, :);
       firstEdgeIndex = size(curIsoEdges, 3) + 1;         
       edgeIndex = firstEdgeIndex;
       curIsoEdges(1, 1, edgeIndex : edgeIndex + 2) = p2;
       edgeIndex = edgeIndex + 3;
       curIsoEdges(1, 1, edgeIndex : edgeIndex + 2) = p3;
       gIsocurvEdges(j, i, firstEdgeIndex : edgeIndex + 2) = curIsoEdges(1, 1, firstEdgeIndex : edgeIndex + 2);
        
    end % FOR nbIsoCurv
end % FOR 1 : nbV


end % Function


function [] = exportIsocurvEdges(vi)
global gIsocurvEdges;

nbIsoCurv = size(gIsocurvEdges, 1);

% Create file to write iso-curves
[fid, Msg] = fopen('isocurves.txt', 'wt');
if ( fid == -1 )
    % something went wrong while opening the file
    error(Msg);
end%if

for i = 1 : nbIsoCurv
    % Write separator
    fprintf(fid, '#\n');
        
    edges = gIsocurvEdges(i, vi, :);
    nbEdges = size(edges, 3) / 6;
    
    for k = 1 : nbEdges        
        p2 = edges(1, 1, k : k + 2);
        p3 = edges(1, 1, (k + 2) : (k + 5) );
        
        % Write next edge of iso-curve to the output file
        fprintf(fid, '[');
        for index = 1 : 3
            fprintf(fid, '%f ', p2(1, 1, index));
        end % FOR
        fprintf(fid, ']');
        fprintf(fid, '[');
        for index = 1 : 3
            fprintf(fid, '%f ', p3(1, 1, index));
        end % FOR
        fprintf(fid, ']\n');
        
    end % FOR nbEdges
    
end % FOR isoCurv

% Close the file
fclose(fid);

end % Function
